library(eplusr)
library(sp)
library(raster)
library(ncdf4)

library(Matrix)
library(energy)
library(FNN)
library(MBC)

month_list <- list('Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec')
root <- "D:/data/GCM_DATA/bias_correction_temp/"
#mod_namelist <- list("MOHC-HadGEM2-ES", "MPI-M-MPI-ESM-MR", "NCC-NorESM1-M")
mod_namelist <- list("MOHC-HadGEM2-ES")
for (month in month_list){
  for (mod_name in mod_namelist) {
    # read calibration period monthly data
    era5_cal_name <- paste(root, "ERA5_cal_period_", month, ".csv", sep = "")
    mod_cal_name <- paste(root, mod_name, "_cal_period_", month, ".csv", sep = "")
    
    era5_cal <- read.csv(era5_cal_name,sep=',',header=FALSE)
    mod_cal <- read.csv(mod_cal_name,sep=',',header=FALSE)
    
    # read validation and projected period monthly data
    mod_val_name <- sub("cal", "val", mod_cal_name)
    mod_val <- read.csv(mod_val_name,sep='',header=FALSE)
    mod_proj_name <- sub("cal", "proj", mod_cal_name)
    mod_proj <- read.csv(mod_proj_name,sep=',',header=FALSE)
    
    # bias correct validation period monthly data using quantile delta mapping method
    fit.qdm_val <- QDM(o.c=era5_cal$V1, m.c=mod_cal$V1, m.p=mod_val$V1, ratio=FALSE)
    qdm.c <- fit.qdm_val$mhat.c
    qdm.p <- fit.qdm_val$mhat.p
    
    write.csv(qdm.c, sub(".csv", "_bc_qdm.csv", mod_cal_name))
    write.csv(qdm.p, sub(".csv", "_bc_qdm.csv", mod_val_name))
    
    # bias correct projected period monthly data using quantile delta mapping method
    fit.qdm_proj <- QDM(o.c=era5_cal$V1, m.c=mod_cal$V1, m.p=mod_proj$V1, ratio=FALSE)
    qdm.pp <- fit.qdm_proj$mhat.p
    
    write.csv(qdm.pp, sub(".csv", "_bc_qdm.csv", mod_proj_name))
  }
}

library(eplusr)
library(sp)
library(raster)
library(ncdf4)

library(Matrix)
library(energy)
library(FNN)
library(MBC)

month_list <- list('Jan', 'Feb', 'Mar', 'Apr', 'May', 'Jun', 'Jul', 'Aug', 'Sep', 'Oct', 'Nov', 'Dec')
root <- "D:/data/bias_correction_temp/"
mod_namelist <- list("MOHC-HadGEM2-ES", "MPI-M-MPI-ESM-MR", "NCC-NorESM1-M")

for (month in month_list){
  for (mod_name in mod_namelist) {
    # read calibration period monthly data
    era5_cal_name <- paste(root, "ERA5_calibration_period_", month, ".csv", sep = "")
    mod_cal_name <- paste(root, mod_name, "_calibration_period_", month, ".csv", sep = "")
    
    era5_cal <- read.csv(era5_cal_name,sep=',',header=FALSE)
    mod_cal <- read.csv(mod_cal_name,sep=',',header=FALSE)
    
    # read validation and projected period monthly data
    mod_val_name <- sub("calibration", "validation", mod_cal_name)
    mod_val <- read.csv(mod_val_name,sep='',header=FALSE)
    mod_proj_name <- sub("calibration", "projected", mod_cal_name)
    mod_proj <- read.csv(mod_proj_name,sep=',',header=FALSE)
    
    # bias correct validation period monthly data using quantile delta mapping method
    fit.qdm_val <- QDM(o.c=era5_cal$V1, m.c=mod_cal$V1, m.p=mod_val$V1, ratio=FALSE)
    qdm.c <- fit.qdm_val$mhat.c
    qdm.p <- fit.qdm_val$mhat.p
    
    write.csv(qdm.c, sub(".csv", "_bc_qdm.csv", mod_cal_name))
    write.csv(qdm.p, sub(".csv", "_bc_qdm.csv", mod_val_name))
    
    # bias correct projected period monthly data using quantile delta mapping method
    fit.qdm_proj <- QDM(o.c=era5_cal$V1, m.c=mod_cal$V1, m.p=mod_proj$V1, ratio=FALSE)
    qdm.pp <- fit.qdm_proj$mhat.p
    
    write.csv(qdm.pp, sub(".csv", "_bc_qdm.csv", mod_proj_name))
  }
}